Mathematical Medicine and Biology: A Journal of the IMA

Compartmental models of infectious disease transmission make assumptions about human behaviors. Specifically, they parameterize interactions across population groups, assumed to have distinct epidemiologically-relevant behavioral patterns, primarily through contact matrices stratified by demographic variables such as age, gender, or socioeconomic status. Although such demographic characteristics are readily measurable, they may inadequately capture the social and psychological forces that govern protective behaviors. Drawing on 20 waves of a national survey conducted throughout the COVID-19 pandemic in the United States, we show that institutional trust - particularly trust in public health agencies, physicians, and hospitals - is a dominant predictor of protective behavior adoption. For mask wearing during periods of strongest pandemic activity, for example, institutional trust explains more behavioral variance across population groups than age, income, education, and partisan affiliation combined. In unadjusted analyses, the difference in protective behavior adoption between individuals with the highest and lowest trust in the CDC was four- to six-fold larger than the corresponding differences by age, income, or educational attainment, and exceeded the difference between Democratic and Republican respondents. This association was institutionally specific (e.g., the relationship attenuates for trust in banks), and behaviorally specific (e.g., trust in the CDC is associated with protective behaviors but not visiting a doctor). The latter suggests that trust modifies voluntary compliance with public health recommendations rather than access to or use of healthcare. We conclude that compartmental models of disease transmission would be substantially improved by incorporating institutional trust as a stratifying variable. We additionally offer a trust-integrated mathematical modeling framework and recommendations for the data infrastructure needed for its implementation.

The osteocyte network in bone is believed to play an important role for how bone tissues sense and respond to mechanical stimulation. Yet, bone adaptation to mechanical loads is often conceptualised as a simple response to mechanical stimuli, such as Wolffs law, which is based on mechanical variables only and takes no account of the cellular basis of mechanosensation. Wolffs law presumes the existence of a reference mechanical stimulus, the mechanical setpoint, above which bone is consolidated, and under which bone is removed. In this paper, we develop a theory of bone tissue sensing and adaptation based on osteocytes to provide new understanding of the role played by osteocyte signals in mechanical adaptation. In this theory, the mechanical setpoint of Frosts mechanostat is explicitly embodied as osteocyte properties involved in mechanotransduction. The mechanical setpoint is allowed to adapt due to the replacement of osteocytes during remodelling, making the setpoint space and time dependent. We propose a mathematical model to implement this new theory of bone adapation and present numerical simulations of this model to explore how mechanobiological response curves (effective Wolffs laws) are modulated by setpoint adaptation during remodelling. By accounting for varying osteocyte populations within bone tissue, we explore bone adaptation under osteocyte disruptions, which is particularly relevant to age-related bone loss. Our model suggests that biological disruptions of remodelling balance cannot always be compensated by mechanical feedback, and that setpoint adaptation during remodelling may have significant observable consequences, such as hysteresis in bone response signatures that resemble lazy zones.

Type 2 diabetes is characterized by progressive aggregation of islet amyloid polypeptide (IAPP) within the islets of Langerhans, a process strongly implicated in beta-cell dysfunction and loss. Although oligomeric IAPP intermediates are widely considered the principal cytotoxic species, the relative contributions of the many biological and kinetic processes governing their formation, clearance, and conversion into fibrils remain poorly quantified. Here, a mathematical model of IAPP aggregation is developed that incorporates the physiology of beta-cell secretion and the microanatomy of the islet, including capillary-mediated clearance, enzymatic degradation, and the kinetics of oligomer and fibril formation within a well-mixed control volume. Building on the hypothesis that oligomers are the major cytotoxic species, the concept of accumulated cytotoxicity is introduced, defined as the time integral of the oligomer concentration, and a systematic sensitivity analysis of this quantity with respect to all model parameters is performed. The results reveal a striking hierarchy: only two parameters, the basal rate of IAPP monomer secretion and the rate constant for spontaneous oligomer dissociation, exert a first-order influence on long-term accumulated cytotoxicity, with dimensionless sensitivities approaching +1 and -1, respectively, while the effect of all other parameters remains subordinate and decays at long times. The model further shows that capillary clearance, owing to the physical exclusion of oligomers from fenestrated capillaries, selectively reduces fibril accumulation and amyloid deposition without affecting oligomer-mediated cytotoxicity, indicating that amyloid area fraction, the standard histological metric of disease severity, may not be a reliable surrogate for cytotoxic burden. The model predicts that approximately 48% of the islet area is replaced by amyloid after 30 years, broadly consistent with histological observations of advanced disease. These findings identify monomer secretion and oligomer dissociation as the most promising therapeutic targets to limit cytotoxic damage in type 2 diabetes and provide a quantitative framework for evaluating candidate intervention strategies.

Cell shape and mechanics are intricately connected and tightly regulated by mechanochemical events including biochemical signaling, cytoskeletal remodeling, and plasma membrane mechanics. While experimental advances in microscopy have shed light on the intricate coordination involved in cell shape change in response to different cues, the ability to conduct three-dimensional simulations in realistic geometries remains an open computational challenge. In this work, we develop a finite-element framework that incorporates advection-diffusion-reaction equations coupled with equations governing the kinematics of a deformable interface representing the cell membrane. We applied this framework to three distinct coupled mechanochemical systems, each governed by geometric partial differential equations, resulting in large deformations of the interface. In all three examples, our simulations revealed the emergence of feedback between cellular signaling, cytoskeletal organization, and cell shape. In our first two sets of simulations, we observed that cell migration and neutrophil protrusion were regulated by membrane tension-mediated feedback. In our final application, we predicted shape changes of a dendritic spine starting from a realistic geometry, and found that the complex shape of the spine gives rise to localized regimes of actin cytoskeleton remodeling not previously observed with idealized geometries. Thus, our finite-element framework allows us to generate new mechanistic insights for biophysical problems.

Microbes exhibit complex dynamic behavior as the result of a large number of biochemical processes, spatial and temporal interactions, environmental variations, and evolutionary pressure. Although significant progress has been achieved in understanding microbial ecological dynamics, multiple open questions remain, including the microscopic mechanisms of growth and the roles of nutrients and stochasticity. In this work, we present a minimal theoretical approach to clarify the link between consumption of resources by microbes and their growth. A stochastic model that accounts for a single microbial species consuming a single type of resource while growing via cell division is studied analytically and via Monte Carlo computer simulations. We identify three distinct dynamical regimes of microbial growth determined by the relative magnitudes of resource uptake and division rates and initial conditions. We also show that stochasticity influences the dynamic behavior when the amounts of microbes or resources are low. The model recovers Monod growth kinetics and provides a mechanistic interpretation of the Monod constant and maximal growth rate. The theoretical framework presented captures a wide spectrum of dynamic behaviors in microbial systems, providing a clearer microscopic picture to explain their underlying complex mechanisms.

Neutrophils move through narrow pores, convoluted channels, and tight spaces in tissue to find infection sites. Their ability to sense weak chemical gradients, undergo directed motion, and solve such path-finding problems rests on internal GTPase signaling circuits that control the front protrusion and rear retraction of the cell. Here we explore several variants of known core polarity circuits, with local and long-ranged negative feedback, including inhibitor downstream of Rac, Rac-Rho antagonism, and effects of membrane tension. The resulting reaction-diffusion (RD) equations for Rac and Rho are then used to simulate protrusion-retractions along the edge of a simulated motile cell. We visualize how cells navigate through narrow tracks with sharp corners and weak chemical gradients in 2D. Our metrics for cell performance include polarity initiation, wall-collision intensity, and track completion. In this way, we expose how Rac and Rho, together with their immediate down and upstream components can fine-tune neutrophil motility through complex environments. Author SummaryWhite blood cells, attracted to sites of infection, migrate through complex tissues to find their target. Such movement requires a balance between robust polarity in one direction versus flexibility in response to spatial cues such as obstacles and sharp turns. Here we use mathematical modeling to explore known intracellular circuits that regulate front protrusion and rear retraction in directed cell migration. We test several such circuits in simulations of cells moving along zigzag tracks with sharp turns. We demonstrate that a basic cell polarity circuit, on its own, has limited success, since cells tend to get trapped in sharp corners. Known modulators of this core, which add local negative feedback, mutual front-back antagonism, and long-range feedback from membrane tension, improve cell performance. A cell with the full front-back-membrane tension regulatory circuit avoids delays due to traps and obstacle collisions, and moves swiftly through a convoluted passage to its target site.

Immune cell engagers (ICE) such as bispecific antibodies (bsAbs), within an immunological synapse, bind and link CD3 on a T cell to a target antigen (TAA) on a cancer cell, forming a trimer (CD3:bsAb:TAA complex). With sufficient trimer numbers within the synapse, the T cell can become activated and promote cancer cell killing. Elranatamab, a CD3-bispecific antibody for multiple myeloma, has received FDA and EMA filing acceptance (August 2023 and December 2023, respectively) adding to a growing list of bsAbs that are treating patients. In the drug development stages of ICE bsAbs, mechanistic modeling approaches are often used to attain a greater quantitative understanding of the modality, preclinically, and provide human pharmacokinetic and efficacious dose predictions to aide in Phase 1 trial design. To date, the majority of ordinary differential equation (ODE) trimer models treat the tumor compartment as well-mixed and trimer formation is governed by a bulk population reaction not accounting for individual synapses. This lack of discrimination can lead to imprecise analysis when analyzing results across E:T ratios using metrics like trimers per T cell or trimers per target cell. To this end we developed an ODE trimer model based on single-synapse complexes (one target cell/one immune cell) with 2D cross-linking trimer formation. We show computationally that the number of trimers per synapse is invariant to the value of the E:T ratio for a given free bsAb concentration, a property that cannot be captured by non-synapse models. A simple demonstration of this discrepancy using the well-known Betts trimer model is presented. We then apply the Betts trimer model coupled to a tumor growth inhibition (TGI) module to show that our synapse-based trimer model is easy to substitute in to model TGI, including the addition of a trimer-per-synapse activation threshold function for cell killing. Overall, our model attempts to balance mechanistic fidelity while limiting the complexity of the model.

Small molecule modifiers whence bound, allosterically, will alter the binding of a macromolecule to one- or more-cognate substrates/partners via conformational and non-conformational changes. Although allostery is inferred directly from empirical data, the mathematical basis of these models, constraints deployed and choice of parameter(s) are not clear. Here, we present and characterize a discrete-to-continuous mathematical model for ensemble distributions of a ligand-interacting macromolecular species across milieux-dependent conformational states and examine its role in the genesis and progression of cooperative binding. The premise, of our model, is a set of occupancy matrices (sparse, binary, strictly delocalized) which can be partitioned by a probability-based hyperparameter into mutually exclusive proper subsets of occupancy matrices with identical multinomial probabilities. Since each subset is canonical with a constituent occupancy matrix, it is characterized by a unique multinomial probability. The inner product of combinatorial pairs of all mutually exclusive subsets of occupancy matrices, with an expression for the summed transitional probabilities (finite differences between unique multinomial probabilities), is the differentiable matrix of strictly positive real-valued numbers for the system of ensemble distributions. Whilst the harmonic mean is presented as a generic solution for a system of ensemble distributions, the row-wise definite integral for each column is the finite union of open intervals (contiguous, strictly monotone) which in tandem with a set of interval-specific and bounded transitional probabilities constitutes a piecewise smooth curve (path-connected-, closed- and compact-set). Our discrete-to-continuous model is phenomenological and able to recapitulate the basic tenets of cooperative binding whilst offering insights into the genesis and progression of the same.

Reaction diffusion (RD) systems play a fundamental role in numerous biochemical and biophysical processes. Here, we present a novel analytical framework for solving RD equations by applying the Wentzel Kramers Brillouin Jeffreys (WKBJ) formalism to Ca nanodomains generated by individual membrane channels, a widely used paradigm for intracellular Ca signaling. Previous models have primarily focused on stationary Ca nanodomains while neglecting diffusion and saturation of intracellular Ca buffers and sensors. In contrast, we derive analytical solutions without these simplifying assumptions. Our analysis demonstrates that sustained Ca influx generates continuously expanding distributions of free Ca, whereas Ca bound buffers and sensors propagate as traveling waves. These predictions are supported experimentally by measurements of one-dimensional fluorescence profiles produced by single-channel activity and two-dimensional profiles generated by whole cell Ca currents. The analytical framework developed here readily extends Michaelis Menten type kinetics to reaction diffusion systems and may therefore be broadly applicable to biochemical and biophysical processes in which diffusion cannot be neglected.

Infectious disease transmission evolves as a dynamic process shaped by biological mechanisms, population behavior, and intervention policies, yet public health responses are often driven by lagging indicators. Accurate short- and long-term disease forecasting is essential for the timely deployment of intervention strategies, healthcare capacity planning, and uncertainty-aware, risk-informed decision-making. To address this challenge, three broad classes of forecasting models have traditionally been used: statistical, machine learning, and mechanistic approaches. However, each of these modeling paradigms faces fundamental limitations. In particular, traditional statistical models often lack the flexibility needed to capture complex disease dynamics, machine learning approaches require large, high-quality data streams, and mechanistic models are notoriously difficult to calibrate. To overcome these challenges, we propose a novel physics-informed machine learning (PIML) framework for forecasting infectious disease dynamics. Our approach simultaneously forecasts new case and hospitalization counts, along with other key epidemiological quantities such as the time-varying reproduction number. This is achieved through the design of a machine learning model and estimation strategy regularized by a system of differential equations that encode disease dynamics of the SIHR model, thereby bridging the gap between purely data-driven and mechanistic models. We demonstrate the proposed methodology through in-depth numerical studies and an application to COVID-19 data collected in the state of South Carolina.

Cancer therapy balances between two competing objectives - treatment efficacy against the tumour and the risk of treatment related severe adverse events, including patient death. Most existing optimal control theory (OCT) formulations rely on optimising heuristic cost functionals that lack direct clinical interpretability. In clinical practice treatment efficacy and patient tolerability are primarily assessed through survival metrics and adverse event rates. Here we introduce the Continuous Lifetime Payoff (CLP), a novel OCT objective functional that directly links treatment decisions to patient survival. It explicitly incorporates tumour dynamics, tumour eradication, and patient mortality from tumour progression, drug-related toxicity and age. We fit age-related mortality from life tables and infer parameters from simulated survival data. The CLP provides a clinically grounded framework for optimising chemotherapy regimens.

Binding reactions in effectively one-dimensional systems, such as proteins diffusing along DNA or other filaments, pose a fundamental coarse-graining challenge because stochastic trajectories are recurrent in one dimension and therefore do not admit a unique, separation-independent macroscopic association rate. As a result, continuum rate equations are not exact in 1D even for initially homogeneous systems. Here we develop a practical framework for mapping stochastic 1D reaction-diffusion dynamics onto effective kinetic models. Using mean-first-passage arguments and particle-based simulations, we define a density-dependent association rate and a corresponding single-rate approximation, and quantify when each provides an accurate description of the underlying stochastic dynamics. We implement 1D reaction-diffusion with excluded volume in the NERDSS software using a free-propagator reweighting algorithm and validate it against known pairwise and many-body limits. Our results show that ordinary rate equations with a single effective rate can accurately reproduce 1D reaction kinetics when the dimensionless parameter governing the ratio of intrinsic to diffusion-limited reactivity is small, with excellent agreement in the strongly rate-limited regime and increasing deviations as diffusion control strengthens. We further show that excluded volume in 1D can appreciably alter both kinetics and equilibrium populations, even at modest particle densities, by reducing accessible length and introducing blockade effects. Together, these results provide quantitative guidance for selecting between spatial simulations, density-dependent rate models, and single-rate continuum descriptions of reversible 1D binding reactions.

The control of uterine activity during pregnancy is a complex process that involves regulating myometrial excitability across multiple scales. While numerous studies have investigated various regulatory mechanisms and established the contributions of ion channels and gap junctions, how these mechanisms interact to produce observed changes in uterine activity remains poorly understood. Pivotal to these efforts are computational models that effectively capture gestational changes in excitability across scales. In this study, we propose a multiscale computational modeling framework that can reproduce measured activity at the cellular and tissue scales at a given gestational stage. At the cellular level, we identify key ion currents underlying the observed electrophysiological properties based on a literature review of their regulation and a sensitivity analysis of the Tong 2011 uterine smooth muscle cell activation model. The conductances of these ion currents are then fit to reproduce characteristic resting membrane potentials and burst properties using Bayesian optimization. To extend to the tissue level, we employ an anisotropic monodomain model, parameterized by the resistivity of late pregnancy uterine muscle, to investigate electrical propagation in a two-dimensional section of uterine tissue. We then apply the multiscale model to study myometrial activation in late pregnancy and elucidate the contributions of ion channel and gap junction regulation in transitioning the uterus from a quiescent state to labor. Our resulting model successfully reproduces measured electrophysiological properties at the cellular level and characteristic single-spike and burst-propagation patterns at the tissue level across the three late-pregnant time points analyzed (days 16/17, 18/19, and 20/21) in a murine model. Furthermore, our results suggest that the regulation of the conductances of the voltage-dependent potassium current (IK1), L-type calcium current (ICaL), and sodium current (INa) is most important in determining preterm uterine excitability. The framework established here will promote the development of more gestationally relevant models to better understand labor progression and the factors involved in dysfunctional labor.

Microglia are the resident immune cells of the central nervous system (CNS), including the brain, spinal cord, and retina, where they serve as the first line of defense against infection and inflammation. Dysregulated microglia activity has been implicated in vision-threatening diseases, highlighting the need to understand how retinal microglia respond to inflammatory stimuli. Importantly, acute inflammation induces substantial redistribution of microglia across retinal layers, yet the mechanisms governing this migration remain poorly understood. Here, we develop the first mathematical model of retinal microglia migration during inflammation to determine how inflammatory exposure, administration route, and species-specific pharmacokinetics shape redistribution dynamics across the retina. The model couples lipopolysaccharide (LPS) pharmacokinetics with microglia migration between the outer plexiform layer (OPL), inner plexiform layer (IPL), and ganglion cell layer/nerve fiber layer (GCL/NFL). Model parameters are calibrated to retinal microglia density measurements from mice following LPS (bacterial endotoxin) challenge, before extending the framework to rats and rhesus macaques to investigate species-specific responses. Simulations also compare how administration route, i.e. intravenous or intraperitoneal injections, alter retinal LPS exposure and subsequent microglia redistribution. Our results suggest that redistribution patterns are driven primarily by LPS delivery route and species-specific pharmacokinetics, rather than the initial microglia distribution across retinal layers. Together, these findings provide new insight into immune cell reorganization in the inflamed retina and demonstrate how mechanistic mathematical modeling can be adapted across experimental designs, administration routes, and animal species.

Vaccination can be a useful intervention for reducing infectious disease burden. Estimating numbers of vaccine-prevented health outcomes is one approach to quantifying the benefits of vaccination. Here we improve a method described by Foppa et al. (1) that assumes vaccination has only direct effects, that is, it cannot prevent infection or onward transmission of the disease. We rederive this method and derive an improved method that increases estimation accuracy with minimal additional analytical complexity. To evaluate the improved method, we simulated disease outbreaks and compared the accuracy of the two methods for estimating prevented disease outcomes. In 84% of simulations performed over a wide parameter space, the improved method had an equal or smaller estimation error compared to the original Foppa method, with 7.9-fold smaller mean error and 44-fold smaller standard deviation of errors. Our study improves a method for estimating prevented burden when assuming vaccination has only direct effects.

Phages must locate and bind to bacterial surface receptors to initiate infection. Their tail fiber configuration critically influences this process. We develop a stochastic model describing surface search as a renewal process, incorporating attachment, detachment, and target-finding steps. Using both numerical simulations and analytical calculations, we quantify how tail fiber number, attachment-detachment rates, and geometric constraints impact the mean and the distribution of time to successful adsorption. Notably, the search efficiency shows a nonmonotonic dependence on tail fibers number, governed by a trade-off between binding stability and diffusion-mediated mobility. This optimum shifts depending on the effective bacterial density, target radius, and fiber reach. Short fiber reach imposes severe geometric constraints, reducing mobility at high tail fiber counts and leading to performance degradation. Our findings suggest that phage adsorption strategies are shaped by a balance between anchoring and exploration, with evolutionary implications for tail fiber design and infection efficiency.

Diffuse gliomas are deadly because the individual tumor cells invade - they travel far from the imageable mass, so it is impossible to remove the tumor completely. On the cellular level, glioma cells seem to be in either a "go" state (in which they do not divide) or a "grow" state (in which they do not migrate). We investigate what this tiny choice has to say about the large-scale speed of the invasion front and whether the implication is sufficiently strong to rule out the classical description of the Fisher-Kolmogorov-Petrovsky-Piskunov (Fisher-KPP) type, in which a single phenotype migrates and proliferates. We derive a two-phenotype reaction-diffusion model with density-dependent switching, and we prove the cooperative (quasi-monotone) structure and the associated comparison principle and study travelling-wave solutions of the model. A leading-edge linearization gives minimal front speed as minimizer of an explicit dispersion relation, and direct simulation verifies the predicted speed. In the experimentally relevant fast switching limit, we find a closed-form expression for the speed, that is, we obtain an effective Fisher-KPP equation with rescaled diffusivity and growth rate, with the fractions of the phenotypes. The "go-or-grow" (GoG) front can move at a maximum speed of half the Fisher speed for the same single-cell motility $D$ and proliferation rate $r$, which occurs only when the cells divide their time equally between the two phenotypes. This bound is directly testable: measurement of the front speed, plus independent determination of $D$ and $r$, discriminates the two hypotheses, and in the GoG case, yields recovery of the phenotype balance. We then extend the result to anisotropic (DTI-informed) invasion along white-matter tracts and discuss implications for understanding clinical measurements of growth rate.

Knowledge of the true number of infections over time is valuable for accurately predicting the future course of an epidemic and planning effective interventions, but the number of cases reported offers only a noisy underestimate of the true number of infections. Disease surveillance strategies based on assessing subsets of the population for current infection (infection prevalence surveys) or antibody presence (seroprevalence surveys) yield crucial information about the number truly infected, but are expensive. To explore impact of survey design considerations--both sample size and sampling frequency--on inference of the number of incident infections over time, we coupled agent-based simulations of respiratory virus epidemics with simulations of infection prevalence and seroprevalence surveys. While returns diminish with increased sample size, we find inference generally improved by increasing survey frequency relative to participants-per-round for any given sample size. After survey rounds reach a sufficient frequency, comparable inference performance may be achieved with either more frequent rounds or more participants per round. Rolling designs with tests conducted each day tend to outperform designs in which testing is divided into discrete rounds. We also show that misspecified assumptions about seroreversion may substantially decrease the quality of inference results.

Insecticide-treated bed nets (ITNs) are central to malaria control. They serve as physical barriers and chemical agents that deter and kill mosquitoes, thereby reducing transmission; however, this form of protection reshapes the immunology of malaria by reducing exposure to Plasmodium parasites and weakening the development of naturally acquired immunity. Against this background, the present study develops a modeling framework to investigate how this tension between protection and immunity alters malaria dynamics once vaccination is introduced as a complementary control strategy and the optimum combination of ITN and vaccine coverage required for malaria control. Unlike standard models that use a fixed proportional reduction in transmission, this study models ITN coverage and efficacy as coupled, time-dependent processes and immunity driven by exposure to infection and vaccination. Rigorous analysis of the model identifies existence conditions for equilibria and shows that malaria can be contained through the synergistic interaction of vaccination, vector control, and immunity-mediated host dynamics. Parameter values of the model are estimated by fitting the model to confirmed malaria case data and the estimated baseline reproduction number using these parameter values is 1.41 (95% confidence interval: 1.34-1.48), confirming sustained transmission. Simulations of the model show that, although ITNs reduce immunity acquisition, their net effect is to reduce infections and improve recovery and survival. Hence, population-level benefits of ITNs outweigh their immunity-reducing effects (particularly when combined with vaccination), leading to a reduced malaria burden. Comprehensive sensitivity analysis indicates that malaria burden is driven mostly by mosquito biting intensity, population capacity, and transmission probabilities; while mosquito mortality, effective treatment, ITN performance, and vaccine efficacy cause the most significant reductions. Additionally, malaria is uncontrollable with universal ITN use and vaccination at baseline efficacy. While individual interventions can achieve control under low transmission, neither 100% ITN coverage nor 100% vaccine coverage can achieve control under moderate-to-high transmission. However, a strong synergy between ITNs and vaccination allows combinations of high efficacy to achieve containment at realistic coverage levels, suggesting that integrated malaria control involving effective vector control, vaccination, and prompt treatment is needed.

Over-activation of mitogen-activated protein kinase (MAPK) signaling underlies numerous inflammatory pathologies that are treated using synthetic glucocorticoids to activate glucocorticoid receptors (GR) and induce expression of dual-specificity phosphatase 1 (DUSP1), which encodes for MAPK phosphatase 1 (MKP1). Despite its importance, the single-cell het-erogeneity of this spatial and temporal pathway has not been fully quantified, several regulatory mechanisms are unclear, and accurate quantitative predictions are not possible based on existing models. To address this challenge, we combined immunocytochemistry (ICC) and single-molecule inexpensive FISH (smiFISH) to quantify endogenous GR transport and DUSP1 transcription dynamics across thousands of single cells following dexamethasone (Dex) stimulation. Using Chemical Master Equations (CME) and likelihood-based inference, we identified clear mechanisms and reaction rates for Dex-driven GR nuclear import; compartment-specific GR degradation; GR-dependent modulation of DUSP1 promoter activation and transcription burst frequencies; DUSP1 transcription, elongation, and transport; and time-dependent and saturation-limited cytoplasmic degradation. Rigorous model comparisons against endogenous, fixed-cell data identifies nuclear GR degradation as the dominant mechanism of receptor clearance, indicates that GR primarily regulates promoter activation, and highlights time-dependent AU-rich element (ARE)-mediated mRNA degradation as a likely mechanism for DUSP1 clearance. With these mechanisms, the fully-parameterized model quantitatively predicts joint distributions of GR translation and decay dynamics, DUSP1 transcription site activity, and nuclear and cytoplasmic mRNA heterogeneity among clonal cells as functions of time and across seven orders of magnitude for Dex induction concentrations. Together, these results show that total DUSP1 mRNA levels emerge from the balance between GR-driven activation and cytoplasmic mRNA decay, with the inferred model quantitatively predicting single-cell distributions across held-out conditions.